{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"../assets/logo.png\" width=\"50\" align=\"left\"> \n",
    "\n",
    "# Doppler\n",
    "\n",
    "***\n",
    "\n",
    "#### Prerequesites\n",
    "- Sampling Data\n",
    "- Range"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Imports\n",
    "import numpy as np                  # Scientific Computing Library\n",
    "import matplotlib.pyplot as plt     # Basic Visualization Library"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Intro - Back to Our Previous Example\n",
    "\n",
    "<img src=\"../assets/speedometer.png\" width=\"400\">\n",
    "\n",
    "Think back to our movie example, we went through some steps to extracting range information from our radar data. This is good for reconstructing the dramatic movie scene, but there are still a lot of things missing. One piece of key information that doesn't necessarily show up on the screen is still left to be discovered. That, of course, is speed. By using radar, we should be able to find out how fast an object is approaching, or fleeting. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## What the Radar Sees\n",
    "\n",
    "Last time, we used the way the ADC samples were captured to find range. We said that since they were captured equispaced apart in time, they formed a time log of sorts we could access. Conceptually, let's think about what would happen if we used that same technique for capturing chirps as well as samples. \n",
    "\n",
    "Some object in front of the radar will reflect portions of the first chirp. Now the second, third, and so on chirps hit this object and reflect. For simplicity, let's assume that the object is stayed within a single range bin for the duration of the entire frame and is moving at a constant velocity. There should be some set of samples across the separate chirps that hold this valuable information. Think about what each of them look like relative to the previous and next."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Doppler Effect\n",
    "\n",
    "The doppler effect is what happens when a wave hits a moving object. The intuition behind it is simple, if you throw a bouncy ball at something moving it will either return faster or slower than normal. For waves, what is affected is the frequency of the wave after it bounces off of the object since the peaks in the wave either spread out or get closer together. This means there is some relation between the frequency of the wave and the velocity of the object it hits (additionally the source of the wave if it isn't already stationary). The elementary equation for this is...\n",
    "\n",
    "$f' = \\frac{v+v_0}{v-v_s}f$\n",
    "\n",
    "However, this equation works with sound waves (which do not obey the equation $c=\\lambda f$ that the EM waves radar generates). We need a new method for finding this change in velocity with waves that are moving at the speed of light."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Discovering Meaning in the Chirps\n",
    "\n",
    "In addition to multiple ADC samples, we have multiple chirps. For an object at range $x$ from the radar, when we receive the respective ADC sample, the product will be a complex number with some phase. If the object is moving away from the radar, the respective ADC sample of the second chirp will come in at a very slightly delayed time. This is because the object also moved slightly away in that miniscule amount of time. Althought this movement is miniscule, the change in phase of the wave can be clearly seen."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Processing Doppler Information\n",
    "\n",
    "***"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 1 - Range FFT\n",
    "\n",
    "Let's start with something that we've already done, but expand it slightly. We know from the last reading that we can do a FFT across the data samples we receive to get range information. However, the data I'm providing this time is from a whole **frame** instead of a single **chirp**. Adjust your method to account for this."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Read in frame data\n",
    "frame = np.load('../assets/simple_frame_1.npy')\n",
    "\n",
    "# Manually cast to signed ints\n",
    "frame.real = frame.real.astype(np.int16)\n",
    "frame.imag = frame.imag.astype(np.int16)\n",
    "\n",
    "# Meta data about the data\n",
    "num_chirps = 128 # Number of chirps in the frame\n",
    "num_samples = 128 # Number of ADC samples per chirp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#################### TODO ####################\n",
    "\"\"\"\n",
    "    Task: Perform the range FFt on the entire frame\n",
    "    Output: range_plot [chirps, range_bins]\n",
    "    \n",
    "    Details: This should be very similar to last time, although the input is now of shape [chirps, adc_samples]\n",
    "\"\"\"\n",
    "\n",
    "# Range FFT\n",
    "\n",
    "##############################################"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Visualize Results\n",
    "plt.imshow(np.abs(range_plot).T)\n",
    "plt.ylabel('Range Bins')\n",
    "plt.title('Interpreting a Single Frame - Range')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 2 - Doppler FFT\n",
    "\n",
    "Looking into the produced plot, it's starting to look like the radar definitely sees something at various ranges. Notably, we can see peak lines at range bins ~40 and ~115. Still, it's hard to tell exactly what. We should now try and also obtain velocity information."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Take a sequential FFT across the chirps\n",
    "range_doppler = np.fft.fft(range_plot, axis=0)\n",
    "\n",
    "# FFT shift the values (explained later)\n",
    "range_doppler = np.fft.fftshift(range_doppler, axes=0)\n",
    "\n",
    "# Visualize the range-doppler plot\n",
    "# plt.imshow(np.log(np.abs(range_doppler).T))\n",
    "plt.imshow(np.abs(range_doppler).T)\n",
    "plt.xlabel('Doppler Bins')\n",
    "plt.ylabel('Range Bins')\n",
    "plt.title('Interpreting a Single Frame - Doppler')\n",
    "plt.show()\n",
    "\n",
    "plt.plot(np.abs(range_doppler))\n",
    "plt.xlabel('Doppler Bins')\n",
    "plt.ylabel('Signal Strength')\n",
    "plt.title('Interpreting a Single Frame - Doppler')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Interesting, there are definitely more things that we can tell from these plots just by decoding the velocity, but we need to know how to decipher them first. This is easy, and very similar to deciphering range from range bins. Additionally, we FFT shifted across the chirps to make it more intuitive to understand. That line in the middle at doppler bin 64 is called zero doppler, meaning everything along that line is static/not moving relative to the radar. This means everything to the left (bins<64) is negative doppler, or moving towards the radar and the opposite for the other half of the doppler bins.\n",
    "\n",
    "\n",
    "Some things that you may have observed:\n",
    "1. Much of the received signal translates to having zero doppler, which makes sense if you think about it because most of the objects around us (and the radar) are not moving and thus zero doppler relative to us.\n",
    "2. The plots show at range bin ~40, there is a grouping of peaks in intensity off to the left, meaning an object is most likely moving towards the radar.\n",
    "3. Also, at range bin ~115, we see there is a peak in the middle of the doppler bins, meaning there is probably a highly reflective static object in front of the radar. These described peaks are more clearly shown in the second plot."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Note\n",
    "\n",
    "It's possible you noticed some things after finishing these two steps. \n",
    "1. Why did we take the range FFT/what if we don't need range?\n",
    "    - You really don't need to, but many tasks do require the range information that comes out of a radar. If yours doesn't need that information then you could probably get away with not taking it.\n",
    "2. Why did we have to take the range FFT before we did the doppler FFT?\n",
    "    - Again, no reason. Futhermore, for this case, the properies of the 2D FFT hold and we can actually do either direction. I can prove it right here..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Range FFT -> Doppler FFT\n",
    "range_bins = np.fft.fft(frame, axis=1)\n",
    "fft_2d = np.fft.fft(range_bins, axis=0)\n",
    "\n",
    "# Doppler FFT -> Range FFT\n",
    "doppler_bins = np.fft.fft(frame, axis=0)\n",
    "rfft_2d = np.fft.fft(doppler_bins, axis=1)\n",
    "\n",
    "print('Max power difference: ', np.abs(fft_2d - rfft_2d).max())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Which is pretty close to zero (not necessarily zero because of bit level rounding errors)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Step 3 - Unit Conversion\n",
    "\n",
    "All the units of the data we produced are of some type of \"bin\". Similarly to range resolution, we have a doppler resolution aka velocity resolution. This can be applied in a similar way, but first let's derive an equation. Again we start with some given equations...\n",
    "\n",
    "- $\\omega = \\frac{4\\pi vT_c}{\\lambda}$ - Rotational frequency of phasor due to object moving at $v$ velocity\n",
    "    - $v$ - Velocity\n",
    "    - $T_c$ - Sampling period\n",
    "    - $\\lambda$ - Wavelength\n",
    "- $\\Delta\\omega \\gt \\frac{2\\pi}{N}$ - Minimum change in rotation of phasor to be resolved by radar\n",
    "    - $N$ - Number of sample points\n",
    "\n",
    "Now its your turn to try and derive our velocity resolution equation. \n",
    "\n",
    "Solve this equation: $\\Delta v > ???$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Data sampling configuration\n",
    "c = 3e8 # Speed of light (m/s)\n",
    "sample_rate = 2500 # Rate at which the radar samples from ADC (ksps - kilosamples per second)\n",
    "freq_slope = 60 # Frequency slope of the chirp (MHz/us)\n",
    "adc_samples = 128 # Number of samples from a single chirp\n",
    "\n",
    "start_freq = 77.4201 # Starting frequency of the chirp (GHz)\n",
    "idle_time = 30 # Time before starting next chirp (us)\n",
    "ramp_end_time = 62 # Time after sending each chirp (us)\n",
    "num_chirps = 128 # Number of chirps per frame\n",
    "num_tx = 2 # Number of transmitters\n",
    "\n",
    "# Range resolution\n",
    "range_res = (c * sample_rate * 1e3) / (2 * freq_slope * 1e12 * adc_samples)\n",
    "print(f'Range Resolution: {range_res} [meters/second]')\n",
    "\n",
    "# Apply the range resolution factor to the range indices\n",
    "ranges = np.arange(adc_samples) * range_res"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Make sure your equation translates to the following\n",
    "velocity_res = c / (2 * start_freq * 1e9 * (idle_time + ramp_end_time) * 1e-6 * num_chirps * num_tx)\n",
    "print(f'Velocity Resolution: {velocity_res} [meters/second]')\n",
    "\n",
    "# Apply the velocity resolution factor to the doppler indicies\n",
    "velocities = np.arange(num_chirps) - (num_chirps // 2)\n",
    "velocities = velocities * velocity_res\n",
    "\n",
    "powers = np.abs(range_doppler)\n",
    "\n",
    "# Plot with units\n",
    "plt.imshow(powers.T, extent=[velocities.min(), velocities.max(), ranges.max(), ranges.min()])\n",
    "plt.xlabel('Velocity (meters per second)')\n",
    "plt.ylabel('Range (meters)')\n",
    "plt.show()\n",
    "\n",
    "plt.plot(velocities, powers)\n",
    "plt.xlabel('Velocity (meters per second)')\n",
    "plt.ylabel('Reflected Power')\n",
    "plt.title('Interpreting a Single Frame - Doppler')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#  Analyze the Data\n",
    "\n",
    "***\n",
    "\n",
    "## Exercise\n",
    "\n",
    "I'm now going to give you some data that has *something* about it that neither range nor velocity themselves can solve. It's your job to try and figure out why this data is different. Doing this task and understanding why it's important will help you see how to make sense of radar data and start to see a bigger picture."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# Read in data\n",
    "frame = np.load('../assets/doppler_example_1.npy')\n",
    "\n",
    "# Manually cast to signed ints\n",
    "frame.real = frame.real.astype(np.int16)\n",
    "frame.imag = frame.imag.astype(np.int16)\n",
    "\n",
    "print(f'Shape of frame: {frame.shape}')\n",
    "\n",
    "# Data configuration\n",
    "num_chirps = 128\n",
    "num_samples = 256\n",
    "\n",
    "num_rx = 4\n",
    "num_tx = 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#################### TODO ####################\n",
    "\"\"\"\n",
    "    Task: Analyze the radar data (range)\n",
    "    Output: range_plot [chirps, range_bins]\n",
    "    \n",
    "    Details: First make some initial observations about the range data.\n",
    "    \n",
    "             Note that this frame actually has 3 different dimensions (axes) to it.\n",
    "             Axis 0 - Chirps\n",
    "             Axis 1 - Virtual Antennas ***ignore the deep meaning of this for now, just treat it as \n",
    "                                          8 copies of [128, 128] frames and work around it then accumulate them together***\n",
    "             Axis 2 - ADC Samples\n",
    "\"\"\"\n",
    "range_plot = None\n",
    "\n",
    "# Range FFT\n",
    "\n",
    "\n",
    "# Take magnitude and then accumulate (along the virtual antennas)\n",
    "\n",
    "\n",
    "##############################################"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot the range plot for each chirp\n",
    "plt.plot(range_plot.T)\n",
    "plt.show()\n",
    "\n",
    "# Feel free to add any additional plots for futher analysis\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#################### TODO ####################\n",
    "\"\"\"\n",
    "    Task: Analyze the radar data (range-doppler)\n",
    "    Output: range_doppler [doppler_bins, range_bins]\n",
    "    \n",
    "    Details: Figure out why this data is special by producing a range-doppler to work with. \n",
    "                (Hint: Look more closely at the peak at range bin ~60)\n",
    "    \n",
    "             Note that this frame actually has 3 different dimensions (axes) to it.\n",
    "             Axis 0 - Chirps\n",
    "             Axis 1 - Virtual Antennas ***ignore the deep meaning of this for now, just treat it as \n",
    "                                          8 copies of [128, 128] frames and work around it then accumulate them together***\n",
    "             Axis 2 - ADC Samples\n",
    "\"\"\"\n",
    "range_doppler = None\n",
    "\n",
    "# Range FFT\n",
    "\n",
    "\n",
    "# Doppler FFT and FFT shift\n",
    "\n",
    "\n",
    "# Take magnitude and then accumulate (along the virtual antennas)\n",
    "\n",
    "\n",
    "##############################################"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot range-doppler plot as an image\n",
    "plt.imshow(range_doppler.T)\n",
    "plt.xlabel('Doppler Bins')\n",
    "plt.ylabel('Range Bins')\n",
    "plt.title('Analyzing the Range-Doppler Heatmap')\n",
    "plt.show()\n",
    "\n",
    "# Feel free to add any additional plots for futher analysis\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The Range-Doppler Plot\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If you correctly analyzed the range-doppler plot, you should have seen that there are multiple objects at the same range. The only reason we are able to distinguish them is because they are moving in opposite directions from one another. Using this information in the movie example, we would could certainly have a lot more information on friend or foe and how quickly we need a decision."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Extras\n",
    "\n",
    "***"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Benefits from mmWave Radars\n",
    "\n",
    "That last exercise may have seemed very straightforward, but there are already a lot of things within your reach. For example, I want to bring up Google Soli, a project that works with radars for gesture recognition ( https://atap.google.com/soli/ ). They use mmWave radars for their high sensitivity to movement and velocity. These properties enable small gestures using a hand or the fingers to be registered as some pattern in the data. Cameras would have a very hard time reading the minute changes from frame to frame, but mmWave radar doesn't have a problem since they are fundamentally different. After training on what specific gestures look like, Soli can classify the gestures as \"change volume\", \"swipe\", and many more. They go over some of the many benefits to doing this on their website and in their demo videos."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Recap\n",
    "\n",
    "This module has gone over how the radar takes advantage of the **doppler effect** to encode **velocity**. Small changes in the phase of ADC samples over multiple chirps can be used to track micro-movements and can then be converted into velocity. Applying this knowledge, we can perform a **doppler FFT** just like the range FFT to get a **range-doppler heatmap**. This heatmap is simple, but can be very powerful at determining a number of things. We showed one use of this by looking at a dopper-range plot with two objects at the same range, but opposite velocities. Range alone would not be able to help since the **two objects would have an additive affect in the range plot**, but the doppler fft gives us a **whole new dimension** to work with."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Looking Forward\n",
    "\n",
    "So far, we have been making a simplifying assumption about the radar we are using. That is because we have assumed that the radar we are using is a simple single receiver (RX) and single transmitter (TX). We can do a lot better. We now have access to radars using **MIMO TX/RX** setups. Just as we've learned that the ADC samples and chirps give us valuable measurements, MIMO gives us a new piece of information, **angle of arrival**."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "***\n",
    "\n",
    "#### Contributors\n",
    "- Dash Kosaka\n",
    "\n",
    "#### Questions, Issues, etc.?\n",
    "Contact by...\n",
    "- email - presenseradar@gmail.com\n",
    "- github - https://github.com/presenseradar/openradar"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
